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DISKS  SURROUNDING  MASSIVE  STARS:  WHEN  COMPUTATIONAL 
MODELS  ARE  CONFRONTED  BY  OBSERVATIONS 

C.  E.  Jones,1  C.  Tycner,2  J.  A.  Benson,3  D.  J.  Hutter,3  H.  R.  Schmitt,4,5  T.  A.  A.  Sigut,1  and  R.  T.  Zavala3 

RESUMEN 

Las  estrellas  masivas  se  encuentran  insertas  en  materia  gaseosa  circunstelar;  a  veces  tambien  se  encuentra 
polvo  y  moleculas.  Aunque  los  discos  son  a  veces  demasiado  pequenos  para  ser  detectados  directamente,  este 
material  puede  ser  detect  ado  en  el  espectro  de  la  radiacion  que  observamos  de  la  estrella.  A  menudo,  el  material 
circunestelar  posee  una  distribution  tipo  disco,  pero  el  proceso  ffsico  que  forma  y  mantiene  estos  discos  no  es 
bien  entendido.  Las  estrellas  Be  son  un  ejemplo  de  rotadores  rapidos,  estrellas  calientes,  cuyo  espectro  en 
longitudes  de  onda  optica  muestra  ambas  lineas  de  hidrogeno  en  emision  y  frecuentemente,  lineas  de  emision 
de  metales  una  vez  ionizados,  clebido  a  la  presencia  del  disco. 

Hemos  calculado  modelos  teoricos  computacionales  de  discos  circunestelares  de  estrellas  Be,  usando  un 
codigo  que  incorpora  un  numero  de  mejoras  sobre  tratamientos  previos  de  la  estructura  termica  del  disco, 
incluyendo  una  composition  quimica  realista.  Estos  modelos  pueden  predecir  los  perfiles  de  lfnea  espectrales 
y  anchos  equivalentes,  la  distribution  espectral  de  energia  y  la  polarization  del  continuo.  Los  modelos  con  es¬ 
tructura  termica  precisa  y  campos  de  radiacion  son  esenciales  para  interpretar  correctamente  las  observaciones. 
Estos  modelos  pueden  tambien  predecir  imagenes  en  el  piano  del  cielo  en  importantes  longitudes  de  onda  y 
por  lo  tanto  estan  idealmente  adecuados  para  compararlos  con  observaciones  interferometricas.  Demuestro  que 
nuestros  modelos  pueden  ser  restringidos  por  comparacion  directa  con  observaciones  interferometricas  opticas 
para  la  region  que  emite  Ha  y  por  perfiles  de  lfnea  Ha  contemporaneos.  Las  comparaciones  detalladas  de  nues- 
tras  predicciones  con  interferometrfa  y  espectroscopfa  de  Ha  entregan  ajustadas  restricciones  a  los  parametros 
libres  del  modelo  para  estos  sistemas  estrella-disco. 

ABSTRACT 

Many  massive  stars  are  embedded  within  gaseous  circumstellar  matter;  sometimes  dust  and  molecules  are  also 
present.  Though  the  disks  are  sometimes  too  small  to  be  detected  directly,  this  material  can  be  detected  in  the 
spectrum  of  radiation  we  observe  from  the  star.  Often,  the  circumstellar  material  has  a  disk-like  distribution, 
but  the  physical  processes  that  form  and  maintain  these  disks  are  not  well  understood.  Be  stars  are  an  example 
of  rapidly  rotating,  hot  stars,  whose  spectra  at  optical  wavelengths  show  both  hydrogen  emission  lines  and, 
frequently,  emission  lines  from  singly  ionized  metals  due  to  the  presence  of  a  disk. 

We  have  computed  theoretical  models  of  circumstellar  disks  for  Be  stars,  using  a  non-LTE  radiative 
transfer  code  which  incorporates  a  number  of  improvements  over  previous  treatments  of  the  disk  thermal 
structure,  including  a  realistic  chemical  composition.  These  models  can  predict  spectral  line  profiles  and 
equivalent  widths,  spectral  energy  distributions,  and  continuum  polarization.  Models  with  accurate  thermal 
structures  and  radiation  fields  are  essential  to  interpreting  observations  correctly.  These  models  can  also  predict 
images  on  the  plane  of  the  sky  in  important  wavelengths  and  are  therefore  ideally  suited  for  comparison  with 
interferometric  observations.  I  will  demonstrate  that  our  models  can  be  constrained  by  direct  comparison 
with  optical  interferometric  observations  for  the  Ha  emitting  region  and  by  contemporaneous  Ha  line  profiles. 
Detailed  comparisons  of  our  predictions  with  Ha  interferometry  and  spectroscopy  place  very  tight  constraints 
on  the  model  free  parameters  for  these  star-disks  systems. 
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1.  INTRODUCTION 

Be  (B-emission)  stars  are  rapidly  rotating  B  stars 
surrounded  by  a  disk- like  distribution  of  gas.  The 
observed  emission  is  due  to  radiative  recombination 
within  the  gaseous  disk.  The  variation  of  the  width 
of  these  emission  lines  as  a  function  of  the  stellar 
rotation  velocity  is  consistent  with  the  assumption 
that  one  is  viewing  an  equatorially  concentrated, 
rapidly  rotating,  slowly  expanding,  disk-like  circum- 
stellar  envelope  at  a  variety  of  inclination  angles  for 
different  stars.  These  stars  also  have  a  high  degree 
of  linear  polarization  due  to  scattering  of  light  off 
ionized  matter,  further  evidence  of  a  non-spherical 
distribution  of  circumstellar  material. 

Despite  decades  of  research,  there  are  no  mod¬ 
els  which  successfully  explain  the  existence  of  this 
circumstellar  material.  Struve  (1931)  proposed  the 
first  model  which  suggested  that  the  presence  of  the 
disk  is  due  to  rapid  stellar  rotation.  This  general 
view  of  a  flattened  circumstellar  disk  has  been  con¬ 
firmed  for  several  stars  by  means  of  direct  imaging  of 
the  circumstellar  disk,  using  either  optical,  infrared, 
and/or  radio  interferometry.  (See  Quirrenbach  et 
al.  1997).  More  recently,  interferometric  observa¬ 
tions  that  can  spatially  resolve  the  disk  are  becoming 
available.  See  for  example  Chesneau  et  al.  (2005); 
McAlister  et  al.  (2005);  Tycner  et  al.  (2006).  We 
know  that  Struve’s  model  is  too  simplistic;  many  ad¬ 
ditional  models  have  been  proposed  which  combine 
rapid  rotation  with  other  mechanisms  such  as  stel¬ 
lar  pulsation  (Rivinius  et  al.  2003),  wind  compressed 
disks  (Bjorkman  &  Cassinelli  2003),  and  magneti¬ 
cally  confined  disks  (Donati  et  al.  2001).  None  of 
these  models  have  been  completely  successful  (Porter 
&  Rivinus  2003). 

Many  of  the  theoretical  models  used  in  the  past 
to  predict  observables  from  these  systems  have,  by 
necessity,  required  the  use  of  a  large  number  of  free 
parameters,  such  as,  disk  opening  angle,  density  dis¬ 
tribution,  temperature,  and  disk  geometry.  For  the 
work  presented  in  this  paper,  the  non-LTE  radia¬ 
tive  transfer  code  of  Sigut  &  Jones  (2007)  is  used  to 
predict  the  flux  in  the  Ha  emitting  region,  the  Ha 
profiles,  and  the  other  predicted  observables.  The 
code  has  a  number  of  improvements  over  previous 
disk  models  available  in  the  literature.  For  example, 
the  disk  has  a  realistic  chemical  composition  that  is 
important  for  computing  heating  and  cooling  rates. 
This  improvement  combined  with  the  enforcement 
of  radiative  equilibrium  results  in  a  self-consistent 
temperature  distribution  for  the  disk.  It  is  crucial 
to  determine  the  state  of  the  gas  accurately,  so  that 
observations  can  be  interpreted  correctly  eventually 


leading  to  the  development  of  successful  dynamical 
models. 

In  this  work,  we  assume  a  power-law  density  dis¬ 
tribution  within  the  disk  in  the  radial  direction,  with 
R  as  the  distance  from  the  centre  of  the  central  star. 
This  density  distribution  can  be  set  by  only  two 
free  parameters,  the  power-law  index,  n,  and  the  as¬ 
sumed  density  at  the  stellar  surface,  p0 ■  The  density 
distribution  in  the  vertical  direction  (perpendicular 
to  the  equatorial  plane),  Z ,  is  approximately  hydro¬ 
static.  This  simple  density  distribution  is  given  by: 

p(R,Z)  =  p0(^y  e-(#)2,  (1) 

where  H  is  the  scale  height  in  the  Z  direction. 
(Please  see  Sigut  &  Jones  2007  for  greater  detail). 

Two-dimensional  images  of  the  predicted  energy 
distribution  in  the  plane  of  the  sky  can  be  produced 
in  any  wavelength  regime  of  interest.  These  images 
can  then  be  transformed  for  direct  comparison  with 
interferometric  observations.  These  comparisons  are 
used  to  place  tight  constraints  on  the  density  param¬ 
eters,  n  and  pa. 

Here,  we  present  a  review  of  our  work  for  the 
Be  stars,  y  Oph  (B2  Ve),  and  ft  Dra  (B6  Illpe). 
The  investigation  for  the  disk  of  \  Oph  resulted  in 
a  best-fit  disk  model  to  be  selected  by  comparison 
with  interferometry  alone.  However,  for  the  case  of 
ft  Dra,  a  preferred  model  could  not  be  selected  on 
the  basis  of  interferometry  alone,  and  contempora¬ 
neous  spectroscopy  in  the  same  wavelength  regime 
is  used  to  remove  the  degeneracy.  This  review  aims 
to  demonstrate,  that  although  the  information  ob¬ 
tained  from  Ha  spectroscopy  and  Ha  interferometry 
is  complimentary,  additional  new  information  can  be 
obtained  by  combining  the  results  of  both  of  these 
techniques. 

2.  OBSERVATIONS 

The  interferometric  observations  were  obtained 
at  the  Navy  Prototype  Optical  Interferometer  in 
Flagstaff  Arizona  in  June  2006  for  %  Oph,  and  in 
February  and  March  2007  for  k  Dra.  We  have  data 
for  multiple  baselines  for  each  star,  ranging  in  length 
from  18.9m  to  53.2m.  The  contemporaneous  spec¬ 
troscopy  was  obtained  at  the  Lowell  Observatory’s 
John  S.  Hall  Telescope.  Please  see  Tycner  et  al. 
(2008)  and  Jones  et  al.  (2008)  for  greater  details  for 
X  Oph,  and  ft  Dra,  respectively. 

3.  RESULTS 

For  each  star,  hundreds  of  models  are  computed 
for  pairs  of  n  and  pQ.  Typically  for  Be  star  disks,  the 
values  of  n  cited  in  the  literature  are  in  the  range 
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2  to  3.5.  This  expected  range  is  based  mainly  on 
previous  results  based  on  models  fitted  to  the  IR 
continuum  (Waters  1986).  In  order  to  ensure  that 
our  search  includes  all  reasonable  pairs  of  n  and  p0, 
we  construct  our  grids  for  n  ranging  from  1.8  to  ~  5 
and  with  pQ  from  ~  10-12  g  cm'3  to  ~  10-8  g  cm'3. 
For  example,  the  grid  of  models  for  k  Dra  can  be 
found  in  Jones  et  al.  (2008),  Figure  1. 

The  integrated  model  intensity  for  each  model  is 
computed  over  a  15nm  spectral  region  centered  at 
Ha  for  consistency  with  the  interferometry.  We  note 
that  this  predicted  intensity  distribution  is  assumed 
to  be  axisymmetric  with  no  over-dense  or  under- 
dense  regions  within  the  disk.  Finally,  the  model  in¬ 
tensities  are  Fourier  transformed  in  preparation  for 
comparison  with  the  interferometry  data.  We  use 
X2  values  to  assess  the  goodness  of  fit  between  the 
models  and  observations. 

3.1.  x  Ophiuchi 

This  bright,  nearby  star  has  pronounced  Ha 
emission  and  is  therefore  well  suited  for  a  study 
of  its  disk  properties  in  this  spectral  region.  The 
comparison  of  ~  500  pairs  of  n  and  pa  with  inter¬ 
ferometry  using  x2  values  for  x  Oph  produced  one 
model  that  was  statistically  better  than  all  of  the 
others.  This  best  fit  model  has  a  reduced  y2  value 
of  1.17  and  corresponds  to  a  model  with  n  =  2.5 
and  p0  =  2.5  x  10“ 11  g  cm-3.  The  predicted  Ha 
line  profile,  corresponding  to  this  model,  is  shown  in 
Figure  1.  We  note  that  the  calculation  of  the  Ha 
profile  and  the  predicted  visibilities  were  the  result 
of  the  same  model  calculation.  The  predicted  line 
profile  was  convolved  with  a  Gaussian  of  FWHM  of 
0.656  A  to  match  the  resolving  power  of  the  obser¬ 
vation. 

This  best  fit  model  is  also  compared  with  ob¬ 
servations  of  the  spectral  energy  distribution  in  the 
infrared  previously  reported  in  the  literature.  See 
Figure  8  of  Tycner  et  al.  (2008).  The  comparison 
shows  reasonable  agreement,  especially  considering 
the  fact  that  x  Oph  is  known  to  be  variable.  Many 
Be  stars  are  known  to  be  variable  on  a  variety  of 
timescales  (Porter  &  Rivinus  2003),  and  this  serves 
to  illustrate  the  need  to  obtain  contemporaneous  ob¬ 
servations. 

3.2.  n  Draconis 

A  similar  procedure,  as  described  above,  was  fol¬ 
lowed  for  k  Dra.  However,  when  the  model  pre¬ 
dictions  were  compared  with  Ha  interferometry  we 
found  that  25  models  had  reduced  y2  <  2.  There¬ 
fore,  further  constraints  were  required  to  select  a  best 


Fig.  1.  Ha  emission  line  for  x  Oph.  The  observation 
was  obtained  on  2006  June  10  (solid  line)  and  has  an 
equivalent  width  of  —7.1  ±  0.2  nm.  The  predicted  line 
(dashed  line)  produced  from  the  best  fit  model  has  an 
equivalent  width  of  —6.6  nm. 

fit  model.  For  all  of  the  models  for  n  Dra,  we  con¬ 
structed  Ha  line  profiles  and  computed  their  corre¬ 
sponding  equivalent  widths.  There  were  14  models 
with  Ha  equivalent  widths  within  2  A  of  the  con¬ 
temporaneous  Ha  observations.  Within  this  subset, 
8  models  had  reduced  x2  <  2  from  comparisons 
with  interferometry.  Therefore,  this  clearly  illus¬ 
trates  that  although  6  models  had  an  Ha  equivalent 
width  that  matched  the  observed  equivalent  width, 
these  did  not  have  the  correct  spatial  distribution 
of  gas  for  agreement  with  interferometry.  The  Ha 
line  profiles  for  the  subset  of  8  were  compared  to  the 
observed  profile  and  we  were  able  to  select  a  best 
fit  model  based  on  both  the  spectroscopy  and  the 
interferometry  comparisons. 

Figure  2  shows  the  observed  and  predicted  Ha 
lines  for  k  Dra  for  the  4  models  that  correspond  to 
the  best  fit  from  both  interferometry  and  Ha  mod¬ 
eling.  Please  see  Jones  et  al.  (2008)  for  a  discussion 
of  the  errors.  We  find  that  the  model  with  n  =  4.2 
and  p0  =  1.5  x  10-10  g  cm-3  is  the  best  fit  in  terms 
of  both  the  spectroscopy  and  interferometry. 

4.  CONCLUSION 

We  have  demonstrated  by  using  detailed  mod¬ 
els  and  interferometric  observations  that  the  spatial 
disk  density  can  be  determined.  Also,  we  showed 
that  the  process  for  finding  the  best-fit  disk  model 
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Fig.  2.  The  observed  Ha  line  for  n  Dra  compared  with 
predicted  lines  for  the  4  models  that  correspond  to  the 
best  fit  from  both  interferometry  and  spectroscopy.  The 
reduced  yft  from  interferometry  for  these  models  are  1.18, 
1.29,  1.18,  and  1.18,  from  left  to  right. 


for  the  Be  stars,  x  Oph  and  k  Dra,  was  necessarily 
different.  For  x  Oph  a  comparison  of  interferometric 
observations  with  models  allowed  a  best  fit  to  be  de¬ 
termined  and  that  this  best  fit  model  matches  other 
observables.  On  the  other  hand,  it  was  necessary  to 
use  both  Ha  interferometry  and  spectroscopy  for  n 
Dra  in  order  to  remove  the  degeneracy  and  find  a 
preferred  model. 

In  future,  we  plan  to  extend  the  work  by  includ¬ 
ing  more  sophisticated  forms  of  the  disk  density  dis¬ 
tribution  eventually  expanding  our  analyzes  to  in¬ 
clude  dynamical  models. 

The  Navy  Prototype  Optical  Interferometer  is  a 
joint  project  of  the  Naval  Research  Laboratory  and 
the  US  Naval  Observatory,  in  cooperation  with  Low¬ 
ell  Observatory,  and  is  funded  by  the  Office  of  Naval 
Research  and  the  Oceanographer  of  the  Navy.  We 
also  thank  the  Lowell  Observatory  for  the  telescope 
time  used  to  obtain  the  Ha  line  spectra  presented 
in  this  review.  This  research  was  supported  in  part 
by  NSERC,  the  Natural  Sciences  and  Engineering 
Research  Council  of  Canada. 

DISCUSSION 

Oliver  Chesneau:  You  have  presented  the  x2 
valley  as  a  function  of  pa  and  n,  and  the  valleys  are 
similar  for  the  Ha  and  the  interferometric  data.  The 
fit  of  the  SED’s  even  though  the  sources  are  variable, 
should  be  considerable  help  the  fitting  process.  The 


determination  of  the  n  factor  is  potentially  important 
for  constraining  the  dynamics  of  the  material. 

Yes,  I  agree  the  value  of  n  is  an  important  diagnostic. 
The  average  value  of  n  could  potentially  be  used  to 
either  support  or  refute  particular  dynamic  models. 

Douglas  Gies:  I  am  impressed  by  the  location  of 
the  cooler  gas  in  your  models:  close  to  the  equatorial 
plane  and  relatively  near  to  the  star.  Do  you  think 
that  this  cooler  gas  might  be  the  source  of  the  shell 
lines  that  are  observed  in  spectra  of  stars  with  incli¬ 
nations  close  to  90°  ?  These  lines  often  correspond 
to  the  spectrum  of  a  star  cooler  than  the  Be  star.  - 
Yes,  the  deep  absorption  troughs  at  the  line  centres 
of  shell  spectra  likely  originate  from  this  cooler  gas 
along  the  line  of  sight. 

D.  Baade  I  would  not  worry  too  much  about  si- 
multaneaously  fitting  the  wings  of  the  Ha  emission 
and  its’  overall  shape.  The  former  are  quite  variable 
and  probably  signal  phases  of  enhavend  star-to-disk 
mass  transfer.  But  a  strong  emission  line  (disk)  of¬ 
ten  shows  little  or  no  response  to  this  at  later  times. 
Therefore,  wings  and  overall  shape  of  the  Ha  emis¬ 
sion  lines  may  not  merit  simultaneous  modeling. 

Chris  Tycner  The  Ha  emssion  was  very  stable 
(at  the  few  percent  level)  during  the  NPOI  interfero¬ 
metric  observations.  No  variations  in  the  wings  were 
detected. 
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